home *** CD-ROM | disk | FTP | other *** search
/ Language/OS - Multiplatform Resource Library / LANGUAGE OS.iso / oper_sys / weyl / weyl_lsp.lha / epolynomial.lisp < prev    next >
Encoding:
Text File  |  1991-11-04  |  7.4 KB  |  228 lines
  1. ;;; -*- Mode:Lisp; Package:Weyli; Base:10; Lowercase:T; Syntax:Common-Lisp -*-
  2. ;;; ===========================================================================
  3. ;;;                  Expanded Polynomials
  4. ;;; ===========================================================================
  5. ;;; (c) Copyright 1989, 1991 Cornell University
  6.  
  7. ;;; $Id: epolynomial.lisp,v 2.9 1991/11/04 16:11:54 rz Exp $
  8.  
  9. (in-package "WEYLI")
  10.  
  11. (defmethod make-epolynomial
  12.     ((domain multivariate-polynomial-ring) compare-function (poly mpolynomial))
  13.   (make-epolynomial domain compare-function (poly-form poly)))
  14.  
  15. (defmethod make-epolynomial
  16.     ((domain multivariate-polynomial-ring) compare-function (poly epolynomial))
  17.   (if (eql compare-function (compare-function poly)) poly
  18.       (make-instance 'epolynomial
  19.              :domain domain
  20.              :compare-function compare-function
  21.              :form (sort (poly-form poly)
  22.                  #'(lambda (a b)
  23.                      (%funcall compare-function
  24.                            (first a) (first b)))))))
  25.  
  26. (defmethod make-epolynomial
  27.     ((domain multivariate-polynomial-ring) compare-function (form list))
  28.   (let* ((dimension (length (ring-variables domain)))
  29.      (evs (get-vector-space (get-lisp-numbers) dimension))
  30.      poly-form)
  31.     (labels ((scan-poly-form (form exp next-var)
  32.            (cond ((poly-coef? form)
  33.               (unless (0? form)
  34.             (loop for i upfrom next-var below dimension
  35.                   do (push 0 exp))
  36.             (push (cons (%apply #'make-element evs (reverse exp))
  37.                     form)
  38.                   poly-form)))
  39.              (t (let ((var-index (poly-order-number form)))
  40.               (loop for i upfrom next-var below var-index
  41.                   do (push 0 exp))
  42.               (setq var-index (1+ var-index))
  43.               (map-over-each-term (poly-terms form) (e c)
  44.                 (scan-poly-form c (cons e exp)
  45.                         var-index)))))))
  46.       (scan-poly-form form () 0)
  47.       (make-instance 'epolynomial
  48.              :domain domain
  49.              :compare-function compare-function
  50.              :form (sort poly-form
  51.                  #'(lambda (a b)
  52.                      (%funcall compare-function
  53.                            (first a) (first b))))))))
  54.     
  55. (defmethod print-object ((poly epolynomial) stream)
  56.   (if (0? poly) (princ 0 stream)
  57.       (with-slots (var-count) poly
  58.         (let ((first? t)
  59.           (variables (ring-variables (domain-of poly))))
  60.       (map-over-each-term (poly-form poly) (e c)
  61.         (cond ((minus? c)
  62.            (setq c (minus c))
  63.            (if (and (1? c) (not (0? e)))
  64.                (princ " -" stream)
  65.                (format stream " - ~S" c)))
  66.           ((null first?)
  67.            (if (and (1? c) (not (0? e)))
  68.                (princ " + " stream)
  69.                (format stream " + ~S" c)))
  70.           (t (unless (and (1? c) (not (0? e)))
  71.                (format stream "~S" c))))
  72.         (loop for var in variables
  73.           for i upfrom 0
  74.           for exp = (ref e i) do
  75.               (unless (lisp:zerop exp)
  76.         (princ " " stream)
  77.         (display var stream)
  78.         (unless (eql 1 exp)
  79.           (format stream "^~D" exp))))
  80.         (setq first? nil))))))
  81.  
  82. (defmethod 0? ((x epolynomial))
  83.   (null (poly-form x)))
  84.  
  85. (defmethod 1? ((x epolynomial))
  86.   (let ((form (poly-form x)))
  87.     (and (null (rest form))
  88.      (0? (le form))
  89.      (1? (lc form)))))
  90.  
  91. ;;; ===========================================================================
  92. ;;;            Polynomial Arithmetic
  93. ;;; ===========================================================================
  94.  
  95. ;; Notice that we can use MAP-OVER-EACH-TERM since it doesn't depend
  96. ;; upon the exponent arithmetic, while UPDATE-TERMS does.
  97.  
  98. (defmacro same-compare-functions ((x y) &body body)
  99.   `(with-slots (compare-function) x
  100.      (unless (eq compare-function (slot-value ,y 'compare-function))
  101.        (error "EPolynomials don't have same compare function: ~S and ~S"
  102.           ,x ,y))
  103.     ,@body))
  104.  
  105. ;; These should really check to make sure that the number of variables
  106. ;; hasn't changed.
  107. (defmethod-binary plus epolynomial (x y)
  108.   (same-compare-functions (x y)
  109.     (bind-domain-context domain
  110.       (make-instance 'epolynomial
  111.         :domain domain
  112.     :compare-function compare-function
  113.     :form (gterms-plus compare-function (poly-form x) (poly-form y))))))
  114.  
  115. (defun gterms-plus (compare-function x y) ;x and y are term lists
  116.   (let ((ans-terms (list nil))
  117.     (terms nil)
  118.     sum)
  119.     (macrolet
  120.     ((collect-term (.e. .c.)
  121.        `(progn (setf (rest terms) (make-terms , .e. , .c.))
  122.          (setf terms (rest terms)))))
  123.       (setq terms ans-terms)
  124.       (loop
  125.        (cond ((terms0? x)
  126.           (cond ((terms0? y) (return (rest ans-terms)))
  127.             (t (collect-term (le y) (lc y))
  128.                (setq y (red y)))))
  129.          ((or (terms0? y)
  130.           (%funcall compare-function (le x) (le y)))
  131.           (collect-term (le x) (lc x))
  132.           (setq x (red x)))
  133.          ((%funcall compare-function (le y) (le x))
  134.           (collect-term (le y) (lc y))
  135.           (setq y (red y)))
  136.          (t (setq sum (+ (lc x) (lc y)))
  137.         (unless (0? sum)
  138.           (collect-term (le x) sum))
  139.         (setq x (red x) y (red y))))))))
  140.  
  141. (defmethod minus ((x epolynomial))
  142.   (let ((domain (domain-of x)))
  143.     (bind-domain-context domain      
  144.       (make-instance 'epolynomial
  145.         :domain domain
  146.     :compare-function (slot-value x 'compare-function)
  147.     :form (gterms-minus (poly-form x))))))
  148.  
  149. (defun gterms-minus (x)
  150.   (map-over-each-term x (e c)
  151.     (collect-term e (- c))))
  152.  
  153. (defmethod-binary difference epolynomial (x y)
  154.   (same-compare-functions (x y)
  155.     (bind-domain-context domain
  156.       (make-instance 'epolynomial
  157.         :domain domain
  158.     :compare-function compare-function
  159.     :form (gterms-difference compare-function (poly-form x) (poly-form y))))))
  160.  
  161. (defun gterms-difference (compare-function x y) ;x and y are term lists
  162.   (let ((ans-terms (list nil))
  163.     (terms nil)
  164.     sum)
  165.     (macrolet
  166.     ((collect-term (.e. .c.)
  167.        `(progn (setf (rest terms) (make-terms , .e. , .c.))
  168.          (setf terms (rest terms)))))
  169.       (setq terms ans-terms)
  170.       (loop
  171.        (cond ((terms0? x)
  172.           (cond ((terms0? y) (return (rest ans-terms)))
  173.             (t (collect-term (le y) (- (lc y)))
  174.                (setq y (red y)))))
  175.          ((or (terms0? y)
  176.           (%funcall compare-function (le x) (le y)))
  177.           (collect-term (le x) (lc x))
  178.           (setq x (red x)))
  179.          ((%funcall compare-function (le y) (le x))
  180.           (collect-term (le y) (- (lc y)))
  181.           (setq y (red y)))
  182.          (t (setq sum (- (lc x) (lc y)))
  183.         (unless (0? sum)
  184.           (collect-term (le x) sum))
  185.         (setq x (red x) y (red y))))))))
  186.  
  187.  
  188. (defmethod-binary times epolynomial (x y)
  189.   (same-compare-functions (x y)
  190.     (bind-domain-context domain
  191.       (make-instance 'epolynomial
  192.         :domain domain
  193.     :compare-function compare-function
  194.     :form (gterms-times compare-function (poly-form x) (poly-form y))))))
  195.  
  196. (defun gterms-mon-times (poly-terms e c)
  197.   (if (0? c) (terms0)
  198.       (map-over-each-term poly-terms (te tc)
  199.     (collect-term (+ e te) (* tc c)))))
  200.  
  201. (defun gterms-times (compare-function x y)
  202.   (let (;; Multiply x by the first term of y.  This is the initial
  203.     ;; term list we will modify.
  204.     (answer (gterms-mon-times x (le y) (lc y))) 
  205.     e c)
  206.     (setq answer (cons nil answer))
  207.     (loop for (e-y . c-y) in (red y)
  208.       for ans = answer do
  209.       (loop for (e-x . c-x) in x do
  210.     (unless (0? (setq c (* c-x c-y)))
  211.       (setq e (+ e-x e-y))
  212.       ;; Find place to insert this term.
  213.       (loop    for red-ans = (red ans) do
  214.         ;; Sure would be nice if the complier recognized and optimized
  215.         ;; the usages of (red ans)
  216.         (cond ((or (terms0? red-ans)
  217.                (%funcall compare-function e (le red-ans)))
  218.            (setf (red ans) (make-terms e c red-ans))
  219.            (return t))    
  220.           ((= e (le red-ans))
  221.            (setf (lc red-ans) (+ (lc red-ans) c))
  222.            (return t))
  223.           (t (setq ans red-ans)))))))
  224.     (loop for ans on answer
  225.       do (when (and (red ans) (0? (lc (red ans))))
  226.            (setf (red ans) (red (red ans)))))
  227.     (red answer)))                 
  228.